Tous les corrigés se trouvent (ou se trouveront) ici : https://curiousml.github.io/
# si pb d'affichage.
#import os
#os.environ["PATH"] += os.pathsep + "C:\\Users\\Kim Antunez\\Anaconda3\\Library\\bin\\graphviz"
#os.environ["PATH"]
import numpy as np
from p_decision_tree.DecisionTree import DecisionTree
import pandas as pd
data_dict = {'Time Match' : ["Morning", "Afternoon", "Night", "Afternoon", "Afternoon", "Afternoon", "Afternoon", "Night", "Afternoon", "Afternoon", "Afternoon", "Afternoon", "Morning", "Afternoon", "Night", "Afternoon"],
'Type Court' : ["Master", "Grand Slam", "Friendly", "Friendly", "Master", "Grand Slam", "Grand Slam", "Master", "Master", "Grand Slam", "Master", "Grand Slam", "Master", "Grand Slam", "Friendly", "Grand Slam"],
'Surface' : ["Grass", "Clay", "Hard", "Mixed", "Clay", "Grass", "Hard", "Mixed", "Grass", "Hard", "Clay", "Hard", "Grass", "Grass", "Hard", "Clay"],
'Best Effort' : ["Yes", "Yes", "No", "No", "Yes", "Yes", "Yes", "Yes", "Yes", "Yes", "Yes", "Yes", "Yes", "Yes", "No", "Yes"],
'Y' : ["Win", "Win", "Win", "Loose", "Loose", "Win", "Win", "Loose", "Win", "Win", "Loose", "Win", "Win", "Loose", "Win", "Win"]}
data_df = pd.DataFrame(data_dict)
data_df
Y = data_df["Y"]
X = data_df[["Time Match", "Type Court", "Surface", "Best Effort"]]
decisionTree = DecisionTree(np.array(X).tolist(),
["Time Match", "Type Court", "Surface", "Best Effort"],
np.array(Y).tolist(),
"gini")
#Here you can pass pruning features (gain_threshold and minimum_samples)
dot = decisionTree.id3(0,0)
#dot = decisionTree.print_visualTree(render=True)
display(dot)
print("System entropy: ", format(decisionTree.entropy))
print("System gini: ", format(decisionTree.gini))
On considère un jeu de données médicales. Chaque exemple correspond à une tumeur du sein. Les variables explicatives portent sur des caractéristiques observées de la tumeur. La variable à prédire indique s’il s’agit d’une tumeur maligne (0, cancer) ou bénigne (1, pas cancer) à partir de caractéristiques géométriques observées. On pourra consulter une description plus détaillée du jeu de données à l’aide de la commande print(data.DESCR).
import numpy as np
import matplotlib.pyplot as plt
from sklearn.datasets import load_breast_cancer
data = load_breast_cancer()
X = data.data
y = data.target
print(data.DESCR)
Question 1. — Séparer le jeu de données en un échantillon d’apprentissage (X_train,y_train) et un échantillon de test (X_test,y_test).
from sklearn.model_selection import train_test_split
X_train, X_test, y_train, y_test = train_test_split(X, y) # test train spit à 75% - 25%
X_train
y_train
Nous allons travailler dans ce TP avec des SVM. Nous pouvons les entraîner de la façon qui suit où l’argument $C$ correspond à l’inverse de l’hyperparamètre de régularisation $\lambda$ vu dans le cours; autrement dit plus C est grand, et moins la régularisation est importante.
Remarque : Différence SVM et perceptron, on veut maximiser la marge.
Question 2. — Sur l’échantillon d’apprentissage, utiliser la validation croisée ainsi que les courbes de validation afin de choisir une valeur pour l’hyperparamètre C. Entraîner à nouveau le prédicteur, qu’on appellera svm, sur tout l’échantillon d’apprentissage avec la valeur de C choisie. Donner son score sur l’échantillon de test. Définir également un array y_test_predict contenant les prédictions du prédicteur obtenu sur l’échantillon de test.
# Pour C = 1
from sklearn.svm import LinearSVC
svm = LinearSVC(C=1).fit(X_train, y_train)
print('Score: ', svm.score(X_test,y_test))
from sklearn.model_selection import validation_curve
from warnings import simplefilter
from sklearn.exceptions import ConvergenceWarning # afin d'ignorer le warning ci-dessus
simplefilter("ignore", category=ConvergenceWarning) # afin d'ignorer le warning de convergence ci-dessus
C_range = np.linspace(0.01,1,50)
#print(C_range)
train_scores, valid_scores = validation_curve(LinearSVC(), X_train, y_train, "C", C_range, cv=10)
train_scores_mean = np.mean(train_scores,axis=1)
valid_scores_mean = np.mean(valid_scores,axis=1)
plt.plot(C_range,train_scores_mean,label="scores d'apprentissage")
plt.plot(C_range,valid_scores_mean,label="scores de validation")
plt.legend()
plt.xlabel('log(C)') #et non C
plt.ylabel('score')
plt.show()
Remarque : le modèle est instable pour C entre 0.01 et 1. Penalisons encore plus !
Le graphique ci-après montre qu'en prenant C relativement petit, le modèle est plus stable et plus performant
C_range = np.logspace(-6,-2,50)
C_range
# C très petit
train_scores, valid_scores = validation_curve(LinearSVC(), X_train, y_train, "C", C_range, cv=10)
train_scores_mean = np.mean(train_scores,axis=1)
valid_scores_mean = np.mean(valid_scores,axis=1)
plt.plot(np.log10(C_range),train_scores_mean,label="scores d'apprentissage")
plt.plot(np.log10(C_range),valid_scores_mean,label="scores de validation")
plt.legend()
plt.xlabel('log(C)') #et non C
plt.ylabel('score')
plt.show()
Après avoir bien paramétré notre hyperparamètre $C$, nous observons une nette améliortation du score !
C_best = C_range[np.argmax(valid_scores_mean)]
svm = LinearSVC(C=C_best).fit(X_train,y_train)
print('Score:', svm.score(X_test,y_test))
y_test_predict = svm.predict(X_test)
pour plus d'information sur le paramètre C : https://stats.stackexchange.com/questions/31066/what-is-the-influence-of-c-in-svms-with-linear-kernel
On peut afficher la matrice de confusion du prédicteur obtenu sur l’échantillon de test comme suit.
à partir de la matrice de confusion, nous pouvons construire plusieurs quantités qui vont nous permettre de mesurer la performance des prédicteurs de façon plus fine qu’avec le simple score. Dans problème présent, il s’agit de détecter si une tumeur est maligne (0) ou bénigne (1).
Il est donc bien plus grave prédire qu’une tumeur maligne est bénigne, que l’inverse. On souhaite construire un prédicteur qui détecte 95% des tumeurs malignes, sans toutefois signaler toutes les tumeurs comme malignes.
Nous introduisons à présent plusieurs notions portant sur les prédictions d’un prédicteur sur un échantillon donné, et qui sont définies par rapport à l’étiquette $y^*$, dite positive, c’est-à-dire celle qu’il est le plus important de détecter (ici $y^*$ = 0 qui correspond à la tumeur maligne).
Un vrai négatif et un faux négatif sont définis de façon similaire, et on note TN et FN leurs nombres respectifs.
Question 3. — Calculer la précision et le rappel à partir des coefficients de la matrice de confusion.
La fonction confusion_matrix prend en argument (dans cet ordre) vrai_y, y_prédit. Le coefficient C[i,j] correspond au nombre d'observations du groupe i prédites comme étant dans le groupe j. Dans le cas de prédiction 0-1 (où 0 est considéré comme négatif, et 1 positif), C[0,0] correspond donc au nombre de vrai négatif, C[0,1] au nombre de faux positif, C[1,0] au nombre de faux négatif, C[0,0] au nombre de vrai positifs. Ces notations correspondent au cadre du test d'hypothèse que vous avez déjà rencontré, ou la réponse est négative si on ne rejete pas H0 ("il ne se passe rien"), et positive si on rejete H0 (par rapport à une hypothèse H1, "il se passe quelque chose") : on détecte qu'un groupe a une moyenne différente de l'autre, que la proportion de cancer est plus importante chez les fumeurs, ou autres exemples auxquels vous pouvez penser. Attention, ce n'est pas le cas ici : la réponse "négative" (hypothèses nulles, il ne se passe rien, la tumeur n'est pas grave) correspond au 1, la réponse positive (à détecter, il se passe quelque chose, la tumeur est maligne) au 0.
from sklearn.metrics import confusion_matrix
confusion_matrix_test = confusion_matrix(y_test,y_test_predict)
print(confusion_matrix_test)
On voudrait minimiser le nombre de tumeurs malignes non détectées (c'est plus grave que l'inverse. )
Rq : Attention en test, un positif correspondrait à Y=1. Mais ici ce n'est pas le cas car les tumeurs malignes qu'on veut détecter sont associées avec une étiquette qui vaut 0. Donc dans la doc de confusion_matrix : elle prend le vrai Y et ensuite le Y prédit. Ici il y a donc une subtilité.
tp = confusion_matrix_test[0,0] ### habituellement c'est les TN c'est que là on
# considère que 0 = positif et 1 = négatif
fp = confusion_matrix_test[1,0]
fn = confusion_matrix_test[0,1]
tn = confusion_matrix_test[1,1]
recall = tp/(tp+fn)
precision = tp/(tp+fp)
taux_faux_positifs = fp/(fp+tn) # proportions de négatifs déclarés positifs
print('Recall :', recall)
print('Precision :', precision)
print('Taux de faux positifs : ',taux_faux_positifs)
Le rappel est plus faible que la précision. Notre classifieur est meilleur pour Y=0 benigne que maligne Y=1 (mais en fait ici maligne est codé 0) (pas terrible). Pas grave d'avoir un taux de faux positif élevé tant que le rappel est bon.
scikit-learn fournit des fonctions qui permettent de calculer directement la précision et le rappel. Cela nous permet de vérifier les résultats obtenus à la question précédente. On notera que l’argument pos_label=0 permet ici de spécifier que l’étiquette considérée comme positive est 0 (tumeur maligne).
from sklearn.metrics import recall_score, precision_score
print('Recall: ',recall_score(y_test,y_test_predict,pos_label=0))
print('Precision: ',precision_score(y_test,y_test_predict,pos_label=0))
print('Accuracy score: ', accuracy_score(y_test, y_test_predict))
Question 4. — Définir une fonction false_positive_rate, qui prend trois arguments: y_true qui contient les étiquettes réelles d’un échantillon, y_predict qui contient les étiquettes prédites par un prédicteur, et pos_label qui spécifie l’étiquette considérée comme positive; et qui renvoie le taux de faux positifs. Calculer le taux de faux positifs donnés sur l’échantillon de test par le prédicteur svm obtenu à la question 2.
Cette fonction calcule le taux de faux positif qui n'existe peut-être pas dans les packages. Il faut remarquer que le taux de faux positif : proportion de négatifs déclarés positifs. négatifs (y_true != pos_pos) et nb de négatifs déclarés positifs (y_predict == pos_label)
def false_positive_rate(y_true,y_predict,pos_label):
return np.sum(y_true[y_predict == pos_label] != pos_label)/np.sum(y_true != pos_label)
print(false_positive_rate(y_test,y_test_predict,0))
On souhaite donc construire un prédicteur dont le rappel est d’au moins 95 %. Notons $(\hat w, \hat b) \in \mathcal{R}^{30}$ les coefficients du prédicteur construit à la question 2. Le prédicteur lui-même s’écrit, puisque les étiquettes sont ici Y = {0,1}
$ f_{\hat w, \hat b} = \begin{cases} 1 & \langle\,\hat w, \hat b\rangle + \hat b > 0 \\ 0 & \langle\,\hat w, \hat b\rangle + \hat b < 0 \end{cases} $
Nous allons ici choisir une famille de prédicteurs de la forme :
$ f_{\hat w, \hat b}^{(\tau)} = \begin{cases} 1 & \langle\,\hat w, \hat b\rangle + \hat b > \tau \\ 0 & \langle\,\hat w, \hat b\rangle + \hat b < \tau \end{cases} $
où $\tau \in R$ est un seuil à choisir. On voit que ces prédicteurs vont avoir tendance, par rapport à $f_{\hat w, \hat b}$, à plus prédire 0 ou plus prédire 1, selon la valeur de $\tau$ .
Question 5. — Définir une fonction modified_predictor qui prend en argument un array contenant les entrées $X_i$ d’un échantillon, ainsi qu’une valeur tau (pour $\tau$ ), et qui renvoie un array contenant les prédictions $f_{\hat w, \hat b}^{(\tau)}(X_i)$. On pourra utiliser la fonction svm.decision_function contenue dans le prédicteur svm, qui prend un argument un array contenant les entrées Xi d’un échantillon, et qui renvoie un array contenant les valeurs $\langle\,\hat w, \hat b\rangle + \hat b$.
def modified_predictor(X,tau):
return (svm.decision_function(X) >= tau).astype('int') # ou *1
Le séparateur détermine un hyperplan tel que :
Mais nos données ne sont pas séparables. Donc on va déplacer l'hyperplan vers la gauche de telle sorte qu'on augmente le rappel (la proportion de positifs non détectés augmente) mais on fait du coup passer des benigne à droite mais moins grave.
On nous demande de modifier le prédicteur pour nous retourner le seuil. On veut fwb de tau où on peut modifier la décision. On garde le même bêta_hat et omega_hat (même hyperplan séparateur) mais on va le déplacer en gardant le même axe. On déplace l'hyperplan du côté des tumeurs benignes pour augmenter notre rappel. On va faire passer des tumeurs malignes de l'autre côté mais c'est pas grave.
Ainsi, si on augmente tau, on va avoir moins de point etiquetées 1 (benigmes) et on va augmenter conjointement :
Pour chaque valeur de $\tau$ , le prédicteur correspondant $f_{\hat w, \hat b}^{(\tau)}$ donne sur l’échantillon d’apprentissage un certain rappel et un certain taux de faux positifs. Pour nous aider à choisir parmi cette famille un prédicteur ayant un rappel d’au moins 95 % et un faible taux de faux positifs, nous allons tracer la courbe ROC : il s’agit de tracer les points de coordonnées (rappel , FPR) (calculés sur l’échantillon d’apprentissage) lorsqu’on fait varier $\tau$ (notons que $\tau$ n’est pas explicitement représenté sur cette courbe).
Question 6. — A quoi ressemble la courbe ROC d’un bon (resp. mauvais) prédicteur?
Lorsqu'on augmente tau, c'est à dire qu'on augmente le nombre de points étiquetés positifs (ou on diminue les négatifs), on s'attend à observer deux phénomènes
Plus précisément, remarquons que rappel/FPR = TP/FP * #{y_true = 1}/#{y_true = 0} Pour un mauvais prédicteur, l'étiquette "P" est décorréllée de la vérité, i.e. du fait que la valeur réelle soit 1 ou 0, et donc le rapport TP/FP est a peu près égal au rapport #{y_true = 0}/#{y_true = 1}. Les points (rappel, FPR) d'un mauvais prédicteur sont donc situés sur la diagonale. Pour un bon prédicteur, le taux de faux négatifs est faible et le rappel est grand. La courbe ROC est donc bien au dessus de la diagnale.
Question 7. — Quelles sont les valeurs minimale et maximale pour $\tau$ qu’il convient de considérer ? Se donner une grille de 100 valeurs à essayer pour $\tau$ . Pour chacune de ces valeurs, calculer sur l’échantillon d’apprentissage le rappel et le taux de faux positifs du prédicteur correspondant $f_{\hat w, \hat b}^{(\tau)}$. Tracer la courbe ROC correspondante.
decision_function_train = svm.decision_function(X_train)
tau_range = np.linspace(np.min(decision_function_train),np.max(decision_function_train),100)
Remarque : Il ne sert à rien de considérer des tau plus petits que le min de svm.decision_function ou plus grand que le max de svm.decision_function sur notre échantillon car tous les points seront du même côté
recalls = []
fprs = []
for tau in tau_range:
y_train_predict = modified_predictor(X_train,tau)
recalls.append(recall_score(y_train,y_train_predict,pos_label=0))
fprs.append(false_positive_rate(y_train,y_train_predict,0))
plt.plot(fprs,recalls)
plt.show()
Très bon prédicteur et au dessus de la diagonale. A chaque point correspond un tau. Quand on augmente tau on se déplace vers la droite (on augmente le nb de positifs donc de faux positifs et donc le rappel)
Question 8. — Parmi les valeurs de $\tau$ considérées, déterminer celle donnant sur l’échantillon d’apprentissage un rappel supérieur à 95 % et le plus faible taux de faux positifs. Pour ce seuil $\tau$, calculer le score, le rappel, la précision et le taux de faux positifs du prédicteur correspondant $f_{\hat w, \hat b}^{(\tau)}$ sur l’échantillon de test.
Objectif :
recalls_array = np.array(recalls)
fprs_array = np.array(fprs)
good_enough_recalls_index = (recalls_array >= .95) #index autorisés : rappel supérieur à 0,95
tau_best = (tau_range[good_enough_recalls_index])[np.argmin(fprs_array[good_enough_recalls_index])]
#on veut index minimum parmi les bons index
print('meilleur tau :', tau_best)
y_test_predict = modified_predictor(X_test,tau_best)
from sklearn.metrics import accuracy_score
print('Accuracy score: ', accuracy_score(y_test, y_test_predict))
print('Recall: ', recall_score(y_test, y_test_predict,pos_label=0))
print('Precision: ', precision_score(y_test, y_test_predict,pos_label=0))
print('False positive rate: ', false_positive_rate(y_test, y_test_predict,pos_label=0))
Précédemment :
l'accuracy est plus faible que le svm précédent mais nous avons un meilleur rappel. Le taux de faux positifs est très élevé
Conclusion : on a entraîné un SVM pour détecter séparation en utilisation validation croisée. Mais on a voulu le modifier. Au lieu d'accepter benigme tous les points dont la fonction de prédiction est négatif, on veut inférieur à tau, et on a choisi le meilleur tau qui vérifiait la contrainte rappel > 0,95. On modifie le classifieur pour augmenter le rappel et ne pas risquer de commettre un certain type d'erreur.
Dans cette 2e partie de TD, on étude des arbres de décision comme dans la partie décision. Cette fois au lieu d'avoir un arbre de décision sur des variables catégorielles, on le fait sur des variables continues.
Pour visualiser les arbres de décision qui vont être construits dans le TP, il convient d’installer les librairies graphviz, python-graphviz et pydotplus. Cela peut être fait en lançant Anaconda Navigator puis en allant sur la page Environments. On considère le jeu de données immobilières. Il s’agit de prédire le prix de biens immobiliers en fonction d’informations concernant l’environnement. Il s’agit donc d’un problème de régression. On pourra consulter la description du jeu de données pour davantage d’informations.
import numpy as np
from sklearn.datasets import load_boston
data = load_boston()
print(data.DESCR)
X = data.data
y = data.target
from sklearn.model_selection import train_test_split
X_train, X_test, y_train, y_test = train_test_split(X,y,random_state=17)
L’argument random_state=17 permet de substituer la valeur 17 à chaque appel d’un générateur de nombre aléatoire lors de l’exécution de la fonction train_test_split. Cela permet de rendre la fonction déterministe et ainsi d’obtenir les mêmes échantillons d’apprentissage et de test à chaque exécution. Cela facilitera la comparaison des différents résultats, car les arbres de décision obtenus dans ce problème dépendent fortement de la composition de l’échantillon d’apprentissage. On peut entraîner un arbre de décision pour ce problème de régression de la façon suivante.
from sklearn.tree import DecisionTreeRegressor
dt = DecisionTreeRegressor(random_state=0).fit(X_train,y_train)
Par défaut, le critère d’impureté employé en régression est l’erreur quadratique moyenne (MSE), et pour chaque noeud, la segmentation est choisie aléatoirement parmi celles minimisant le critère. Pour chaque arbre de décision entraîné dans ce TP, nous allons cependant imposer un comportement déterministe à l’aide de l’argument randomstate=0. On peut accéder à la profondeur et au nombre de noeuds de l’arbre dt via dt.tree.maxdepth et dt.tree.node_count respectivement.
Question 9. — Définir une fonction tree_summary qui prend un argument un arbre de décision entraîné et qui affiche son score d’apprentissage, son score de test, sa profondeur et son nombre de noeuds. Appliquer cette fonction à l’arbre dt construit ci-dessus. Commenter les résultats. On peut limiter la profondeur de l’arbre à l’aide de l’argument max_depth.
def tree_summary(dt):
#help(dt) pour savoir les sous listes de l'objet.
print("Score d'apprentissage : ",dt.score(X_train,y_train))
print("Score de test : ",dt.score(X_test,y_test))
print("Profondeur de l'arbre : ",dt.tree_.max_depth) #arbre du classifieur et on récupère le descripteur max_depth
print("Nombre de noeuds : ",dt.tree_.node_count) #arbre du classifieur et on récup node_count.
tree_summary(dt)
Le score d'apprentissage vaut 1. Ce n'est pas etonnant car cela arrive souvent avec des variables continues. On peut toujours séparer nos points quitte à avoir des arbres avec une feuille à la fin, on peut décrire parfaitement les données d'apprentissage.
Il est bien plus élevé que le score de test : il y a donc un fort sur-apprentissage. En limitant la profondeur de l'arbre, on arrête de diviser nos cases en sous-cases. On s'assure d'avoir un certain nombre de points donc on diminue la variance et le sur-apprentissage (mais on risque d'augmenter le biais).
=> On va choisir le meilleur trade-off biais variance (limiter sur et sous apprentissage).
Question 10. — Entraîner un arbre de décision, nommé dt2, dont la profondeur est limitée à 3. Observer ses scores et ses caractéristiques. Commenter.
Exemple de sous-apprentissage.
dt2 = DecisionTreeRegressor(random_state=0,max_depth=3).fit(X_train,y_train)
tree_summary(dt2)
Arbre de profondeur 3. Le nombre de noeuds est bien moindre que dans l'arbre précédent.
Le score d'apprentissage a diminué (puisqu'on a augmenté le biais), et le score de test a diminué également ! On a réduit trop drastiquement le nombre de degrés de liberté du modèle, et on observe un phénomène de sous-apprentissage. On peut choisir la profondeur optimale par validation croisée.
On peut visualiser l’arbre de décision dt2 de la façon suivante.
Quelles sont les variables explicatives qui interviennent dans l’arbre de décision dt2 ?
from sklearn.externals.six import StringIO
from IPython.display import Image
from sklearn.tree import export_graphviz
from pydotplus import graph_from_dot_data
foo = StringIO()
export_graphviz(dt2,out_file=foo,impurity=False)
graph = graph_from_dot_data(foo.getvalue())
Image(graph.create_png())
Les variables explicatives intervenant dans l'arbre sont X_12, X_5, X_7, X_10 et X_0, qui correspondent respectivement à LSTAT, RM, DIS, PTRATIO et CRIM (faire un summary pour avoir des détails).
Dans la feuille en bas 3e à droite, il y a encore pleins de données qu'on pourrait encore séparer.
Question 12. — à l’aide d’une 10-validation croisée, choisir le meilleur paramètre pour la profondeur maximale de l’arbre de décision.
from sklearn.model_selection import GridSearchCV
param_grid = {'max_depth':[4,5,6,7,8]}
predictor= GridSearchCV(DecisionTreeRegressor(random_state=0),cv=10,param_grid=param_grid)
predictor.fit(X_train,y_train)
print('Paramètre sélectionné:',predictor.best_params_)
print('Score d\'apprentissage: ',predictor.score(X_train,y_train))
print('Score de test: ',predictor.score(X_test,y_test))
from sklearn.model_selection import validation_curve
depth_range = np.arange(1,17)
train_scores, valid_scores = validation_curve(DecisionTreeRegressor(), X_train, y_train, "max_depth", depth_range, cv=10)
train_scores_mean = np.mean(train_scores,axis=1)
valid_scores_mean = np.mean(valid_scores,axis=1)
plt.plot(depth_range,train_scores_mean,label="scores d'apprentissage")
plt.plot(depth_range,valid_scores_mean,label="scores de validation")
plt.legend()
plt.xlabel('profondeur')
plt.ylabel('score')
plt.show()
La profondeur va de 1 à 17 (on avait 17 dans un arbre non contraint donc on aurait pas d'arbre plus profond). On peut utiliser le critère du rasoir pour choisir la profondeur maximal (par exemple égale à 6 ici). Entre deux modèle équivalent, tu prends le moins complexe.
Profondeur = profondeur maximal (sans surprise) et 81 noeuds.
On voit que le score de validation augmente fortement, est maximal pour maxdepth=6, puis diminue. On entraine donc un arbre dt3 avec une profondeur maximale de 6
dt3 = DecisionTreeRegressor(random_state=0, max_depth = 6).fit(X_train,y_train)
tree_summary(dt3)
Les questions suivantes ont pour but la construction d’un prédicteur obtenu par agrégation de plusieurs arbres de décisions entraînés en faisant varier l’échantillon d’apprentissage.
Plus précisément, en notant $S = (X_i, Y_i)_{i \in[m]})$ l’échantillon d’apprentissage initial et n le nombre d’arbres de décision qu’on souhaite agréger, pour chaque $k \in [n]$, on tire un échantillon de taille m: $S^{(k)} = (X_i^{(k)}, Y_i^{(k)})_{i \in[m]})$ où chaque exemple $(X_i^{(k)}, Y_i^{(k)})$ est tiré uniformément et indépendamment parmi les exemples de S (il s’agit donc d’un tirage avec remise); on entraîne ensuite avec $S^{(k)}$ un arbre de décision $\hat f^{(k)}$ et le prédicteur final $\hat f$ est obtenu en moyennant :
$\forall x \in X \quad \hat f(x) = \frac{1}{n} \sum_{k=1}^n \hat f^{(k)}(x)$
Question 13. — Définir une fonction r2_score qui prend un argument deux arrays y_true et y_predict contenant respectivement, pour un échantillon donné, les vraies étiquettes et les étiquettes prédites par un prédicteur, et qui revoie le score R2 correspondant.
def r2_score(y_true,y_predict):
return 1-np.sum((y_true-y_predict)**2)/np.sum((y_true-np.mean(y_true))**2)
Question 14. — à l’aide de la fonction resample, tirer (avec remise) à partir de l’échantillon d’apprentissage (X_train, y_train) un autre échantillon (Xtrain, ytrain) de même taille (le tirage étant avec remise, certains exemples peuvent apparaître plusieurs fois dans le nouvel échantillon, d’autres vont en être absents). Entraîner à l’aide de ce nouvel échantillon un arbre de décision (avec random_state=0 et les autres paramètres par défaut) et afficher son score R2 en utilisant la fonction r2_score (on vérifiera qu’elle donne le même résultat que la fonction .score intégrée au prédicteur).
from sklearn.utils import resample
X_train_,y_train_ = resample(X_train,y_train)
dt_ = DecisionTreeRegressor(random_state=0, max_depth = 6).fit(X_train_,y_train_) #on peut enlever max_depth ?
tree_summary(dt_)
print('Score de test:',r2_score(y_test,dt_.predict(X_test)))
print(dt_.score(X_test, y_test))
Question 15. — Construire 5 prédicteurs de la même façon qu’à la question précédente en ré-échantillonnant (Xtrain, ytrain) à chaque fois. Calculer les prédictions des prédicteurs obtenus sur l’échantillon de test.
dts_predictions = np.zeros((5,y_test.size))
for n in range(0,5):
X_train_,y_train_ = resample(X_train,y_train)
dts_predictions[n,:] = DecisionTreeRegressor(random_state=0).fit(X_train_,y_train_).predict(X_test)
print('Score de test:',r2_score(y_test,dts_predictions[n,:]))
print()
J'entraîne 5 arbres avec un échantillon d'entraînement différent à chaque fois (resample). Score de l'arbre varie fortement jusqu'à 1/10e. Le score dépend donc de l'échantillon d'entraînement, on a tout intérêt à agréger nos arbres pour diminuer la variance.
Question 16. — On considère le prédicteur défini comme la moyenne des 5 prédicteurs construits à la question précédente. Calculer ses prédictions sur l’échantillon de test et calculer le score R2 correspondant.
dts_aggregated_predictions = np.mean(dts_predictions,axis=0)
print("Score R² du prédicteur agrégé:",r2_score(y_test,dts_aggregated_predictions))
Je calcule la moyenne de la prédiction. Non seulement le score des arbres agrégés est meilleur que le score de l'arbre 0 et l'arbre 3, mais aussi meilleur que le score des meilleurs arbres.
Quand on fait des forêts aléatoires, on réduit non seulement la variance mais on peut aussi obtenir des scores meilleurs que le meilleur des scores.
Par contre, l'inconvénient de cette méthode... L'avantage des arbres de décisions c'est qu'ils sont très visuels. On peut utiliser ça pour expliquer comment pourquoi dans certains quartiers le loyer est plus cher que dans d'autres. Avec l'agrégation, on a de meilleurs prédictions mais moins facilement interprétables.
Question 17. — Définir une fonction tree_aggregation qui prend en argument un entier n_trees, qui entraîne n_trees arbres de décision comme précédemment, et qui affiche le score R2 du prédicteur agrégé sur l’échantillon de test. Exécuter cette fonction avec différentes valeurs pour n_trees et tenter d’obtenir un meilleur score. La procédure qui consiste à construire un nouvel échantillon d’apprentissage à l’aide de tirages uniformes avec remise s’appelle le bootstrap. Celle qui consiste à agréger plusieurs prédicteurs obtenus par bootstrap s’appelle le bagging (pour bootstrap aggregating). Enfin, dans le cas particulier des arbres de décision, le bagging est également appelé forêt aléatoire. Scikit-learn propose un algorithme de forêt aléatoire prêt-à-l’emploi qu’on peut utiliser comme suit.
def tree_aggregation(n_trees):
dts_predictions = np.zeros((n_trees,y_test.size))
for n in range(0,n_trees):
X_train_,y_train_ = resample(X_train,y_train)
dts_predictions[n,:] = DecisionTreeRegressor(random_state=0).fit(X_train_,y_train_).predict(X_test)
dts_aggregated_predictions = np.mean(dts_predictions,axis=0)
print("Score de test:",r2_score(y_test,dts_aggregated_predictions))
tree_aggregation(10)
tree_aggregation(100)
tree_aggregation(1000)
Question 18. — Comparer les résultats obtenus à l’aide de la fonction tree_aggregation et ceux de l’algorithme fourni par scikit-learn.
from sklearn.ensemble import RandomForestRegressor
RandomForestRegressor(n_estimators=1000).fit(X_train,y_train).score(X_test,y_test)
Notre procédure et celle scikit-learn donnent des résultats très similaires à l'aléa intrinsèque près de notre fonction.
Conclusion : on a entraîné un arbre de décisions de (profondeur 17, 709 noeuds) => Surapprentissage. Dans une première partie, on a voulu limiter la profondeur maximale. On a choisi la profondeur maximale par validation croisée. Quand on a choisi la meilleure profondeur, on a pu calculer notre prédicteur. On a pu montrer ensuite que le score de ce prédicteur dépend beaucoup du jeu d'entraînement. On peut le voir en modifiant le jeu d'entraînement. Pour obtenir un prédicteur moins variable, on a entraîné différents arbres et agrégé les prédictions. On obtient un prédicteur à variance réduite non seulement mais aussi meilleur. Au début, on a agrégé 5 prédicteurs et en fait on peut le faire avec un nombre arbitraire. Les random forest diminuent la variance. n_tree est un paramètre à choisir, le plus grand possible (dépend du temps à allouer). donc si on augmente suffisamment n_tree on peut enlever max_depth =6 car très bon. L'inconvénient est qu'on perd en interprétation.